3
;;;;;;;;;;;; BEGIN INTERNAL METHODS OF SCI-NUM

(define sci-num-marker #(sci-num))

;;; EXPORT
(define (sci-num? n)
  (and (pair? n) (eq? (car n) sci-num-marker)))

(define (unsafe-make-scientific coef exponent)
  ;; This function will not be exported, because it allows to creation
  ;; of invalid and imbalanced sci-nums.
  (list sci-num-marker coef exponent))

;;; EXPORT
(define (make-scientific coef exponent)
  ;; Allowed range of exponent is -100000000000000 to +100000000000000
  (cond
   ((or (nan? exponent)
        (nan? coef)
        (and (infinite? exponent)
             (positive? exponent)
             (zero? coef))
        (and (infinite? coef)
             (infinite? exponent)
             (negative? exponent)))
    (unsafe-make-scientific +nan.0 +inf.0))

   ((or (zero? coef)
        (< exponent -100000000000000))
    (unsafe-make-scientific 0.0 -inf.0))

   ((or (infinite? coef)
        (> exponent 100000000000000))
    (unsafe-make-scientific (* +inf.0 coef) +inf.0))

   (else
    (let ((floor-exponent (exact (floor exponent)))
          (fractional-remains (- exponent (floor exponent))))
      (rebalance-sci-num
       (if (= fractional-remains 0)
           (unsafe-make-scientific (inexact coef) floor-exponent)
           (unsafe-make-scientific
            (* (inexact coef) (expt 10 fractional-remains))
            floor-exponent)))))))

(define (rebalance-sci-num n)
  ;; Helper for addition. Converts an invalid sci-num to a valid one
  ;; (e.g. 105e3 to 1.05e5)
  (let loop ((coef (abs (sci-coef n)))
             (exponent (sci-exponent n)))
    (cond
     ((zero? coef)
      (unsafe-make-scientific 0.0 -inf.0))
     ((nan? coef)
      (unsafe-make-scientific +nan.0 +inf.0))
     ((and (infinite? coef) (negative? coef))
      (unsafe-make-scientific -inf.0 +inf.0))
     ((and (infinite? coef) (positive? coef))
      (unsafe-make-scientific +inf.0 +inf.0))
     ((< coef 1)
      (loop (* coef 10)
            (- exponent 1)))
     ((>= coef 10)
      (loop (/ coef 10)
            (+ exponent 1)))
     (else
      (unsafe-make-scientific
       (if (sci-negative? n)
           (- coef)
           coef)
       exponent)))))

;; Conversions

;;; EXPORT
(define (number->sci-num n)
  (make-scientific n 0))

;;; EXPORT
(define (sci-num->number n)
  (* (sci-coef n)
     (expt 10 (sci-exponent n))))

(define (coerce-nan n)
  (cond
   ((sci-num? n) n)
   ((number? n)
    (if (or (nan? n)
            (not (real? n)))
        +nan.0
        n))
   (else
    (error "Expecting a number; got:" n))))

(define (coerce-sci-num n)
  ;; if n is a number, convert it to a sci-num
  ;; if n is a sci-num, keep it
  ;; otherwise, raise error
  (cond
   ((sci-num? n) n)
   ((number? n) (number->sci-num (coerce-nan n)))
   (else
    (error "Expecting a number; got:" n))))

;;; EXPORT
(define (sci-coef n)
  (list-ref (coerce-sci-num n) 1))

;;; EXPORT
(define (sci-exponent n)
  (list-ref (coerce-sci-num n) 2))

;;;;;;;;;;;; END INTERNAL METHODS OF SCI-NUM

;;; EXPORT
(define (sci-zero? n)
  (zero? (sci-coef n)))

;;; EXPORT
(define (sci-infinite? n)
  (and (infinite? (sci-exponent n))
       (not (zero? (sci-coef n)))
       (positive? (sci-exponent n))))

;;; EXPORT
(define (sci-nan? n)
  (or (nan? (sci-coef n))
      (nan? (sci-exponent n))))

;;; EXPORT
(define (sci-negative? n)
  (negative? (sci-coef n)))

;;; EXPORT
(define (sci-positive? n)
  (positive? (sci-coef n)))

;;;; Addition and Subtraction

(define (match-sci-num-exponent n target)
  ;; Helper for addition. Brings the exponent of a sci-num up to another.
  ;; This function assumes that target is greater than the starting exponent
  (let loop ((coef (abs (sci-coef n)))
             (exponent (sci-exponent n)))
    (if (< exponent target)
        (loop (/ coef 10)
              (+ exponent 1))
        (unsafe-make-scientific
         (if (sci-negative? n)
             (- coef)
             coef)
         exponent))))

(define (naive-sci-add a b)
  ;; Add two sci-nums that were previously matched
  (make-scientific
   (+ (sci-coef a) (sci-coef b))
   (sci-exponent a)))

(define (sci-add-two a b)
  ;; Add two sci-nums
  (define-values (less more)
    (let ((c1 (sci-coef a))
          (c2 (sci-coef b))
          (e1 (sci-exponent a))
          (e2 (sci-exponent b)))
      (cond
       ((< e1 e2)
        (values a b))
       ((> e1 e2)
        (values b a))
       ((< c1 c2)
        (values a b))
       (else
        (values b a)))))

  (cond
   ;; First we must deal with edge cases.
   ((or (sci-nan? less)
        (sci-nan? more)
        (and (sci-infinite? less)
             (sci-infinite? more)
             (or (and (sci-positive? less)
                      (sci-negative? more))
                 (and (sci-negative? less)
                      (sci-positive? more)))))
    +nan.0)
   ((or (sci-zero? less)
        (sci-infinite? more))
    more)
   (else
    (naive-sci-add more (match-sci-num-exponent less (sci-exponent more))))))

;;; EXPORT
(define (sci+ . args)
  (let loop ((result 0)
             (next args))
    (if (null? next)
        result
        (loop (sci-add-two result (car next))
              (cdr next)))))

;;; EXPORT
(define (sci- first . args)
  (if (null? args)
      (make-scientific
       (- (sci-coef first))
       (sci-exponent first))
      (let loop ((result first)
                 (next args))
        (if (null? next)
            result
            (loop (sci-add-two result (sci- (car next)))
                  (cdr next))))))

;;;;;;; Multiplication and Division

(define (sci-mult-two a b)
  (make-scientific (* (sci-coef a)
                      (sci-coef b))
                   (+ (sci-exponent a)
                      (sci-exponent b))))

;;; EXPORT
(define (sci* . args)
  (let loop ((result 1)
             (next args))
    (if (null? next)
        result
        (loop (sci-mult-two result (car next))
              (cdr next)))))

;;; EXPORT
(define (sci/ first . args)
  (if (null? args)
       (make-scientific
        (/ (sci-coef first))
        (- (sci-exponent first)))
      (let loop ((result first)
                 (next args))
        (if (null? next)
            result
            (loop (sci-mult-two result (sci/ (car next)))
                  (cdr next))))))

;;;;;; Other Conversions

(define (log10 n)
  (/ (log n) (log 10)))

(define (log10->sci-num n)
  ;; Allows representation of large sci nums
  (cond
   ((and (number? n)
         (not (real? n)))
    (make-scientific +nan.0 +inf.0))
   ((infinite? n)
    (if (positive? n)
        +inf.0
        0))
   ((nan? n)
    +nan.0)
   (else
    (let ()
      (define exponent (exact (floor n)))
      (define coef (expt 10 (- n exponent)))
      (make-scientific coef exponent)))))

(define (sci-log10 n)
  ;; Allows representation of large sci nums
  (cond
   ((sci-zero? n)
    -inf.0)
   ((sci-negative? n)
    +nan.0)
   (else
    (+ (sci-exponent n) (log10 (sci-coef n))))))

(define (split-string-e s)
  ;; Helper for string->sci-num. Splits an exponential number
  ;; by its E or e, and filters trailing and leading whitespace

  ;;; STATES:
  ;; start, sign, integral, point, fractional, e, exponent, esign, trailing, halt

  (define reader-state 'start)
  (define result #f)
  (define in (append '(#\space) (string->list s) '(#\space)))
  (define (reader-next!)
    (case reader-state
      ((start)
       (case (car in)
         ((#\space #\newline #\tab) #f)
         ((#\+ #\-)
          (set! result (list (car in)))
          (set! reader-state 'sign))
         ((#\.)
          (set! result (list '(#\0) #\+))
          (set! reader-state 'point))
         ((#\0 #\1 #\2 #\3 #\4 #\5 #\6 #\7 #\8 #\9)
          (set! result (list (list (car in)) #\+))
          (set! reader-state 'integral))
         (else
          (set! reader-state 'halt))))
      ((sign)
       (case (car in)
         ((#\.)
          (set! result (cons (list #\0) result))
          (set! reader-state 'point))
         ((#\0 #\1 #\2 #\3 #\4 #\5 #\6 #\7 #\8 #\9)
          (set! result (cons (list (car in)) result))
          (set! reader-state 'integral))
         (else
          (set! result #f)
          (set! reader-state 'halt))))
      ((integral)
       (case (car in)
         ((#\.)
          (set! reader-state 'point))
         ((#\0 #\1 #\2 #\3 #\4 #\5 #\6 #\7 #\8 #\9)
          (set! result (cons (cons (car in) (car result))
                             (cdr result))))
         ((#\e #\E)
          (set! result (cons (list #\0) result))
          (set! reader-state 'e))
         ((#\space #\newline #\tab)
          (set! result (cons '(#\0) (cons #\+ (cons '(#\0) result))))
          (set! reader-state 'trailing))
         (else
          (set! result #f)
          (set! reader-state 'halt))))
      ((point)
       (case (car in)
         ((#\0 #\1 #\2 #\3 #\4 #\5 #\6 #\7 #\8 #\9)
          (set! result (cons (list (car in)) result))
          (set! reader-state 'fractional))
         (else
          (set! result #f)
          (set! reader-state 'halt))))
      ((fractional)
       (case (car in)
         ((#\0 #\1 #\2 #\3 #\4 #\5 #\6 #\7 #\8 #\9)
          (set! result (cons (cons (car in) (car result))
                             (cdr result))))
         ((#\e #\E)
          (set! reader-state 'e))
         ((#\space #\newline #\tab)
          (set! result (cons '(#\0) (cons #\+ result)))
          (set! reader-state 'trailing))
         (else
          (set! result #f)
          (set! reader-state 'halt))))
      ((e)
       (case (car in)
         ((#\0 #\1 #\2 #\3 #\4 #\5 #\6 #\7 #\8 #\9)
          (set! result (cons (list (car in)) (cons #\+ result)))
          (set! reader-state 'exponent))
         ((#\- #\+)
          (set! result (cons (car in) result))
          (set! reader-state 'esign))
         (else
          (set! result #f)
          (set! reader-state 'halt))))
      ((esign)
       (case (car in)
         ((#\0 #\1 #\2 #\3 #\4 #\5 #\6 #\7 #\8 #\9)
          (set! result (cons (list (car in)) result))
          (set! reader-state 'exponent))
         (else
          (set! result #f)
          (set! reader-state 'halt))))
      ((exponent)
       (case (car in)
         ((#\0 #\1 #\2 #\3 #\4 #\5 #\6 #\7 #\8 #\9)
          (set! result (cons (cons (car in) (car result))
                             (cdr result))))
         ((#\space #\newline #\tab)
          (set! reader-state 'trailing))
         (else
          (set! result #f)
          (set! reader-state 'halt))))
      ((trailing)
       (case (car in)
         ((#\space #\newline #\tab) #f))))
    (set! in (cdr in)))

  (let loop ()
    (if (or (eq? reader-state 'halt)
            (null? in))
        (if result (reverse result) result)
        (begin (reader-next!) (loop)))))

;;; EXPORT
(define (string->sci-num s)
  (define parsed-value (split-string-e s))
  (if parsed-value
      (make-scientific
       (string->number
        (string-append
         (string (list-ref parsed-value 0))
         (list->string (reverse (list-ref parsed-value 1)))
         "."
         (list->string (reverse (list-ref parsed-value 2)))))
       (string->number
        (string-append
         (string (list-ref parsed-value 3))
         (list->string (reverse (list-ref parsed-value 4))))))
      #f))

;;; EXPORT
(define (sci-num->string n)
  (string-append
   (number->string (sci-coef n))
   "e"
   (number->string (sci-exponent n))))

;;;;;;;;; Other mathematical functions

(define (sci-infinite= a b)
  ;; helper for sci=, tests whether a and b are infinite and equal
  (and (sci-infinite? a)
       (sci-infinite? b)
       (or (and (sci-positive? a)
                (sci-positive? b))
           (and (sci-negative? a)
                (sci-negative? b)))))

;;; EXPORT
(define (sci-expt b p)
  (cond
   ((sci-zero? p)
    1)
   ((and (= (sci-coef b) 1)
         (= (sci-exponent b) 0))
    1)
   ((sci-negative? b)
    (cond
     ((and (integer? p)
           (even? p))
      (sci-expt (sci- b) p))
     ((and (integer? p)
           (odd? p))
      (sci- (sci-expt (sci- b) p)))
     ((and (sci-infinite? b)
           (finite? p)
           (negative? p))
      0.0)
     (else
      +nan.0)))
   (else
    (log10->sci-num (* p (sci-log10 b))))))

;;; EXPORT
(define (sci-sqrt n)
  (sci-expt n 1/2))

;;; EXPORT
(define (sci-log n)
  (/ (sci-log10 n) 0.43429448190325176))

;;; EXPORT
(define (sci-exp n)
  (sci-expt 2.718281828459045 n))

;;; EXPORT
(define (sci< a b)
  (< (sci-coef (sci- a b)) 0))

;;; EXPORT
(define (sci<= a b)
  (or (sci-infinite= a b)
      (<= (sci-coef (sci- a b)) 0)))

;;; EXPORT
(define (sci= a b)
  (or (sci-infinite= a b)
      (= (sci-coef (sci- a b)) 0)))

;;; EXPORT
(define (sci> a b)
  (> (sci-coef (sci- a b)) 0))

;;; EXPORT
(define (sci>= a b)
  (or (sci-infinite= a b)
      (>= (sci-coef (sci- a b)) 0)))

